Role of interactions in time-of-flight expansion of atomic clouds from optical lattices 
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We calculate the effect of interactions on the expansion of ultracold atoms from a single site of an 
optical lattice. We use these results to predict how interactions influence the interference pattern 
observed in a time of flight experiment. We find that for typical interaction strengths their influence 
is negligable, yet that they reduce visibility near a scattering resonance. 
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I. INTRODUCTION 

One of the most important probes of cold atom sys- 
tems is time-of-flight imaging. Turning off all trapping 
potentials, a cloud of cold atoms expands for tens of mil- 
liseconds, and an absorption image is taken. In the far 
field limit, the resulting image can be directly interpreted 
as the momentum distribution of the original cloud, if in- 
teractions among the atoms can be neglected during the 
expansion. Here we critically evaluate the validity of ne- 
glecting such interactions during the expansion from an 
optical lattice. 

The question of how to interpret time-of-flight images 
is crucial. These images have been used, for example, 
to distinguish the Mott insulating and superfluid phases 
[1, 2]. They have also been proposed as a tool to de- 
tect vortices in rotating condensates [3] , and are a crucial 
component of more sophisticated probes such as modu- 
lation spectroscopy [4] and Bragg/Raman spectroscopy 
[5]. 

Interactions between cold neutral atoms are parame- 
terized by the s-wave scattering length a, which is typi- 
cally on the order of 5 — 15nm [1,6]. The scattering length 
is almost always very small compared to the distance A/2 
between sites in an optical lattice, A/2 ps 426nm, where 
A is the wavelength of the laser used to create the opti- 
cal lattice [1]. The scattering length can, however, be- 
gin to approach the size a r of the atomic states in one 
well. For example er r w 75nm for a modest optical lattice 
with depth V ~ 10E R , where E R = S, 2 (27r) 2 /(2mA 2 ) is 
the recoil energy of the lattice. Thus when a few par- 
ticles occupy a single site, their interactions are signifi- 
cant [7] . Experiments have measured the resulting energy 
shifts [8] , and recently used them to study atom number 
statistics [9]. While these on-site interactions are impor- 
tant, by the time the wave-packets have expanded enough 
to overlap with neighboring sites, interactions are greatly 
attenuated. 

Hence in our analysis we include interactions between 
atoms expanding from the same site, but neglect all inter- 
site interactions. Thus we are able to investigate whether 
interactions during the initial expansion period affect the 
interference image. Our estimate of the role of interac- 
tions is a lower bound; there may be further interaction 



effects during later stages of the expansion. Most impor- 
tantly, interference effects could lead to strongly inter- 
acting high density regions at intermediate times [10]. 

Within our approximation, the density profile of the 
many-body system during time of flight depends only on 
the t = wavefunction, and the time dependence of a 
cluster of particles expanding from a single site. In the 
next Sec. II we consider the expansion from a single site. 
Time-of-flight interferometry is considered in Sec. III. As 
a numerical example we consider the expansion of a two- 
dimensional, harmonically trapped cloud forming a su- 
perfluid in Sec. IV. We summarize our results in Sec. V. 

In 2008, Gerbier et al. [11] reported the results of a 
very similar calculation, however they gave very few de- 
tails. More recently, Fang, Lee, and Wang [12] reported 
a complementary investigation, where they used a trun- 
cated Wigner approximation to investigate the role of 
interactions during time-of-flight expansion. Restricting 
the expansion to one dimension (ID), they considered 
the dynamics of 10 atoms released from a 10-site optical 
lattice. As we discuss in section II C, interactions play a 
much larger role in ID expansion than in 3D, and Fang et 
al. consequently found nearly a factor of two attenuation 
of the central Bragg peak compared to the noninteracting 
gas. Using very similar parameters, we find that interac- 
tions during 3D expansion only lead to a 5% reduction 
in the amplitude of the central Bragg peak. 



II. SINGLE SITE EXPANSION 



A. Statics 

An optical lattice is typically modeled as a potential 
of the form 

V(x,y,z) = Vo[sin 2 (£;:E) + sin 2 (fcy) + sin 2 (fcz)] (1) 

where k = 2n/X. Near the local minima one may 
approximate the sinusoidal as a harmonic potential 
V e g = mw 2 r 2 /2, with small oscillation frequency oj r — 
(Tik 2 1 'm) -y/VuT^R- The single-particle ground state in 
this potential is a Gaussian 
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where of = h/muj r . We model the interaction among the 
particles as a contact interaction, 



■ J d 3 rft(r)ft(r)ip(r)ip(r), (3) 



where g = 4-Kh 2 a/m. 

We are concerned about how these interactions modify 
the few-body wavefunction on a single site, and how this 
influences the time-of-flight expansion. As long as a < 
a r , the effects of interactions are captured by a Gaussian 
variational ansatz 

$Ar({r Q }) =e l?JV </> cm (r cm ) J| 4> N {r a - rp) (4) 
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where we have introduced the center of mass r cm = 
and generalized to arbitrary spatial di- 



mension d. For N = 1, Eq. (4) reduces to the harmonic 
oscillator wavefunction 0i (r) . We will use this wavefunc- 
tion as a time-dependent variational ansatz to describe 
the dynamics, hence we have introduced the parameters 
/3o and (3m which arc nonzero only if the cluster is ex- 
panding or shrinking. 

For convenience we will work with dimensionless quan- 
tities, using units where U = m = A/2 = 1. 



B. Dynamics 

One produces a variational estimate of the dynamics 
by minimizing the action 



S = dt ]Jd c 
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We use the trial wave- function in Eq. (4), for which the 
spatial integrations can be performed analytically, and 
allow all variational parameters to be arbitrary functions 
of time. A similar approach has been used to describe 
the role of interactions in the dynamics of a harmon- 
ically trapped BEC, where the atom number is much 
larger [13]. Minimizing the action leads to a second order 
differential equation for the width ctjv, 
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The center-of-mass width cr obeys the same equation, 
but with 5 = 0. Note that in the noninteracting limit all 



N dependence drops out. The sole contribution from the 
optical lattice is the term a^/af. During time-of-flight 
expansion, the optical lattice as well as the harmonic 
trapping potential are removed, and this term no longer 
appears in the equations of motion. 

At time t = we set (3 = 0, and take a to be given by 
the static solution with dt<J = 0. Analytic solutions to 
the resulting algebraic equation can only be found when 
d = 2, where 
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The center-of-mass width is simply o-q = 07 in all dimen- 
sions. 

In two dimensions we can analytically integrate the 
equations of motion, 



(8) 



Ng 



a N (t,d = 2) = (I+2J) *o(t), (9) 



1/4 



A>(t) = MM = 2) 



1 t 



(10) 



2 af. + t 2 ' 

The expressions for cto and (3 apply in all dimensions. 

C. Phase accumulation 

The phase of the expanding cluster is crucial for de- 
termining the observed interference pattern. In terms of 
cttv, one finds 
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The contribution in square brackets arises from the in- 
teractions. Interactions increase the width of the initial 
state, which reduces the contribution of the kinetic en- 
ergy and retards the phase relative to the noninteracting 
expansion. This should be contrasted with the contri- 
bution from the potential energy as well as the direct 
interparticle interactions themselves, which increase the 
energy and advance the phase. To determine the net 
sign of the interaction correction is thus not straightfor- 
ward. In particular, a^ is generally larger than <7o and 
the quantity in square brackets does not have a definite 
sign. 

We produce a rough estimate of the phase accumu- 
lated by replacing a^ with cto in this expression. The 
interaction contribution to the phase will then scale as 
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Whereas the contribution is logarithmically divergent 
in the one-dimensional case, it very quickly reaches a 
finite value in the two-dimensional as well as the three- 
dimensional case. This indicates that the influence of 
interactions is confined to the very beginning of the time- 
of-flight expansion, t < of, essentially corresponding to 
the time required for the cluster to expand to less than 
twice its initial size, where <r(t) 2 w 1/2. (Recall, we are 
using units where lengths are measured in terms of the 
lattice spacing and times, up to numerical constants, in 
units of the inverse recoil energy.) Typically, this means 
that for d = 2 or d = 3 interactions become irrelevant 
well before the clusters overlap. Conversely, interactions 
between clusters can not be neglected during d = 1 ex- 
pansion. 

The non-interacting contribution to the phase is 
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In two dimensions, where we have analytic expressions 
for apf(t), we further find 
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The fact that interactions only modify the prefactor is a 
reflection of the scaling symmetry of the expanding cloud 
in d = 2. 



III. TIME-OF-FLIGHT IMAGES 

Having calculated the expansion dynamics of a single 
cluster of particles, we now explore the consequences for 
the atom density seen in a time-of-flight expansion exper- 
iment. Neglecting correlations between sites, we assume 
that the initial state can be written as a generalization 
of the standard Gutzwiller Ansatz, 

Ms I OO - \ 
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where i runs over all Af s lattice sites Rj , n is the number 
of particles on a given site, and 

|{r Q }> = -l=^(r 1 )^(r 2 )..4Hr n )\0). (16) 



The state is normalized when the norm of the /-vector 
is one, 
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The n-particle wavefunction on site i, $„({r Q , — Rj}), is 
given by Eq. (4). 



Within our approximation, where we neglect interac- 
tions between atoms on different sites, the time evolution 
of Eq. (15) simply amounts to separately time evolving 
each cluster, as described in Sec. II. The resulting density 
profile is 

(*l^ (r )^( r )|*) = ^;^; n |/. n |2|0=( r _R,)|2 (i 8 ) 
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with 



K(r)| 2 = J n^ 2 d 3 r a |$„(r,r 2 ,...,r„)| 2 (19) 
0n,n-i(r)= / K^rMr,^...)^^...). 

In the noninteracting limit both contractions reduce to 
the noninteracting single particle wavefunction, so that 
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In this case one can write a more readily interpretable 
expression for the density [3], 

<*ni|n(r)|*»i) =|0i (r)| 2 [(N N c ) + |A(r)| 2 ] (21) 



(22) 



Here N is the total number of particles in the lattice. 
1 4>i (r) | 2 is a simple gaussian with width ao(t)/V2. Cor- 
rections to the featureless gaussian peak, 7V C and A(r), 
signal the presence of superfluid order in the system. 
is the condensed number of particles, whereas Ofj is the 
expectation value of the annihilation operator on site i 
and thus the superfluid order parameter in the system. 



(23) 



(24) 



Gerbier et at [11] have pointed out that in Eq. (22) above, 
for experimentally relevant expansion times on the order 
of tens of milliseconds, it is necessary to keep the Fresnel 
like terms quadratic in Hi. In the absence of the Fresnel 
terms, the shape of Bragg peaks is simply the Fourier 
transform of the superfluid order parameter. 

Here we go beyond the approximations in Eqs. (20)- 
(24), and include the effects of interactions on the ex- 
pansion. These interactions have two effects. First they 
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broaden each of the expanding clusters. This broadens 
the incoherent background, but it also reduces the con- 
trast of the Bragg peaks. This latter effect occurs because 
of the reduced overlap between the expanding clusters 
with different numbers of particles. Second, the inter- 
actions introduce a nonlinear phase difference between 
the different particle number clusters. This dephasing 
further reduces the contrast of the Bragg peaks. 

The broadening of the incoherent background is quan- 
tified by 
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Clearly a^.cS > Co? reflecting the larger size of the inter- 
acting cluster. 

The influence of interactions on the Bragg peaks is 
quantified by identifying the difference from the nonin- 
teracting wavefunction </>i(r), 

= Mr)e lSiN \SN N \e X p[r 2 Ss N }. (27) 

The overall phase S£n, the width (<Jq 2 — Ss)^ 1 ' 2 , as well 
as the corresponding prefactor \SNn\, affect the peaks. 
The expressions for each of these terms are complicated, 
with 
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In the nonintcracting case the only phases contribut- 
ing to the interference come from the terms fi,nf* n -i 
[cf. Eq. (23)]. Here there are additional contributions as 
given by Eq. (28). 



IV. NUMERICAL EXAMPLE 

To illustrate our results, we consider a two dimensional 
harmonically trapped gas of 87 Rb in optical lattices with 
V = 10.7, 8, 5E R , yielding v r w 20.3, 17.5, 13.9kHz. For 
each of these lattice depths we use a different scattering 
length, taking a = 5, 15.5, 50nm, to give the same ratio 




FIG. 1: (Color online) Initial evolution of the widths of dif- 
ferent number states N = 1 through 10 for a — 50nm and 
Vo = 5-Br. The smallest widths are reached for N = 1 (lowest 
curve), the largest for iV = 10 (highest curve). Interactions 
cause states with larger particle number to broaden and ex- 
pand faster. Times are measured in units of to = m(A/2) 2 /Ti, 
which is approximately 0.25ms for 87 Rb. The range plotted 
is much shorter than a typical time-of- flight experiment. All 
subsequent expansion is ballistic. 



U /At « 15.5, well on the superfluid site of the superfluid- 
Mott transition [14]. Adjusting the chemical potential /io 
in the center of the trap to obtain the same total number 
of particles thus yields identical initial states. To find the 
initial fi. n 's of Eq. (15) we solve the discrete variational 
Gutzwiller problem, minimizing ('Jdii/rJ'i'G), with 



H L =- t^2(a}a,j + a'jai) + Un^hj - 1) 



(31) 
(32) 



where denotes nearest neighbor sites, annihi- 

lates a boson at site i, hi = afa^, and t and U are 
extracted from the non-interacting Wannier wavefunc- 
tions [14]. The corrections to U from using the many- 
body wavefunctions on each site are very small at this 
lattice depth, and the corrections to t are at most 10% 
[7] . We take the expansion to be three dimensional, treat- 
ing the individual wells as spherically symmetric. 

We produce initial conditions by solving Eq. (6) with 
the conditions dta = and /3 = 0. Starting from these 
initial conditions, we numerically integrate Eqs.(6) and 
(11). We then plot the densities, Eq. (18). 

Figure 1 shows the time evolution of the widths er/v of 
the clusters expanding from sites with different particle 
numbers. As one can see, and as discussed in Sec. II, the 
expansion very quickly becomes ballistic. 

Figure 2 shows the time evolution of the phase differ- 
ences — <5£at, see Eq. (28). One sees that when a = 50nm 
the phase difference between clusters of different particle 
numbers are on the order of 27r, and hence the interfer- 
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FIG. 2: (Color online) Overall phase difference between 
wavefunctions with particle occupation differing by one [see 
Eq.(28)], plotted for a = 50nm and Vo = 5Er. Interactions do 
not affect the lowest (blue) curve, corresponding to n = {0, 1}, 
but progressively affect the higher pairs n — {1, 2}, {9, 10}, 
yielding the largest effect for n = {9, 10}, the top (violet) 
curve. The effect is approximately linear in TV. 



ence pattern will be influenced by the interactions. For 
typical 87 Rb parameters, a — 5nm, the phase difference 
is correspondingly smaller. 

Figure 3 shows cuts through illustrative density images 
along the lattice direction after a 100ms time-of-flight. 
This time was chosen to minimize distortions from Fres- 
nel terms [11]. Weak interactions, a = 5nm, have neg- 
ligable effect on the image. While stronger interactions 
a = 50nm begin to reduce the amplitudes of the interfer- 
ence peaks, the peaks remain clearly visible. 

Comparing expansion images at fixed U /At and fixed 
time of flight results in an interesting structure. In order 
to have the same U/At, the initial wavefunctions in the 
stronger interacting case must be larger, resulting in a 
slower initial expansion. Therefore in Fig. 3 we see that 
the central peak is larger for stronger interactions, while 
the satellite peak is smaller. Using our intuition from the 
noninteracting expansion, one can think of this effect as 
being due to the envelope of the Bragg peaks which falls 
off on a scale inversely proportional to the size of the ini- 
tial Wannier states. In the inset of Fig. 3 we normalize 
out this effect by multiplying with the inverse envelope, 
(7r<7Q(t)) 3 / 2 exp[r 2 /cro(i)]. Taking this normalization into 
account, interactions reduce the amplitude of the inter- 
ference peak. The reduction is about 5% for a = 5nm, 
rising to 33% for a = 50nm. As shown in the inset of 
Fig. 3, the reduction is greater for the central Bragg 
peak than for the first satellite peak. 

There are a number of ways of increasing the impor- 
tance of the interactions during time of flight. For ex- 
ample, changing the geometry of the lattice sites influ- 
ences how the cloud expands and how long interactions 
remain relevant: the expansion from needle shaped sites 
is predominantly in the x-y plane, and the 2D scaling 
in Eq. (12) approximately holds. Additionally we have 
studied what happens when one suddenly increases the 
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FIG. 3: (Color online) Density of a nonrotating atomic cloud 
expanded for 100ms, with distance d measured from its center 
along a lattice direction. The interference peaks are clearly 
visible and remain so even in the strongest interacting limit 
considered here. Plotted are a = Onm (blue), 5nm (violet), 
15.5nm (brown), and 50nm (green). These progress from top 
to bottom in the inset and in the satellite peak at larger d, but 
bottom to top in the central peak. The inset shows the cen- 
tral and the first satellite peak with adjustment for different 
expansion velocities, as described in the main text. 



scattering length while releasing the atoms from the op- 
tical lattice. This allows one to independently control 
the size of the initial Wannier states and the scattering 
length. In the expansion shown in Fig. 3, where t/U was 
fixed while changing a, the initial Wannier states were 
larger when a was made larger. 

Starting from the equilibrium state with a — 5nm 
and Vq — 10.7 En, we investigate the expansion for 
a = 15.5nm and a = 50nm. Although we do not show the 
results here, we find that the suppression of the central 
interference peak is roughly a factor of 1.5 greater than 
what is seen in the inset of Fig. 3. 



V. SUMMARY AND DISCUSSION 

We have considered the effect of two-particle interac- 
tions on the time of flight images of cold atoms on optical 
lattices. We show that on-site interactions can be impor- 
tant for these images, but argue that one can neglect the 
interactions between atoms on different sites. 

We find that even if interactions are increased by a 
factor of ten from their normal strength, no qualitatively 
new features appear in the time of flight images. How- 
ever, the quantitative size of the peaks is sensitive to 
the interactions. Given the wide tunability achievable 
by employing Feshbach resonances [15], it is conceivable 
that experiments can study the role of interactions dur- 
ing time of flight. The analysis presented here will fail 
when the scattering length a becomes comparable to the 
size of the Wannier state oy . 

Conceptually it is worth noting that in many electronic 
mesoscopic systems a situation markedly different from 
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the one here is encountered. There the dynamics is de- 
termined by impurity scattering or scattering off system 
boundaries, and at low temperatures are not affected by 
the interaction. Consequently there exists a regime in 
which interactions effectively only add additional phases 
to the relevant propagation amplitudes. Such a regime 
is not identifiable in the system we considered here. In- 
stead, we find that whenever the interactions produce 
relevant phases, they also perturb the dynamics. 
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